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Abstract 

We consider numerical approximations of stochastic differential equations by 
the Euler method. In the case where the SDE is elliptic or hypoelliptic, we show 
a weak backward error analysis result in the sense that the generator associated 
with the numerical solution coincides with the solution of a modified Kolmogorov 
equation up to high order terms with respect to the stepsize. This implies that 
every invariant measure of the numerical scheme is close to a modified invariant 
measure obtained by asymptotic expansion. Moreover, we prove that, up to 
negligible terms, the dynamic associated with the Euler scheme is exponentially 
mixing. 

Resume 

Nous etudions la discretisation d'une equation differentielle stochastique (EDS) 
par le schema d'Euler. Dans le cas d'une EDS elliptique ou hypoelliptique nous 
montrons un resultat d'analyse d'erreur retrograde : une fonctionnelle de la solu- 
tion numerique est proche de la solution d'une equation de Kolmogorov modifiee 
a des ordres arbitrairement eleves par rapport au pas de discretisation. On ob- 
tient ainsi que toute mesure invariante du schema numerique est proche d'une 
mesure invariante modifiee obtenue par developpement asymptotique. De plus, 
le schema est exponentiellement melangeant a des ordres arbitrairement eleves. 

Keywords: backward error analysis, stochastic differential equations, exponen- 
tial mixing, numerical scheme, Kolmogorov equation, weak error. 
MSG number: 65C30, 60H35, 37M25 
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1 Introduction 



In the last decades, backward error analysis has become one of the most powerful 
tool to analyze the long time behavior of numerical schemes applied to evolution 
equations. The main idea can be described as follows: Let us consider an ordinary 
differential equation of the form 

m = fivit)), 

where / : M" is a smooth vector field, and denote by (p{{y) the asso- 

ciated flow. By definition, a numerical method defines for a small time step 
r an approximation of the exact flow (pr'- We have for bounded y € R", 
^r(y) = ^t{y) + 0{t'^'^^) where r is the order of the method. 

The idea of backward error analysis is to show that can be interpreted as 
the exact flow ipt'' or a modified vector field defined as a series in powers of r 

/r = / + T^r + T^+Vr+l + • • • , 

where f£,i>r are vector fields depending on the numerical method. In general, 
the series defining does not converge, but it can be shown that for bounded y, 
we have for arbitrary N 

<l>r{y) = (y) + Cnt'' , 

where is the truncated series: 

Under some analyticity assumptions, the constant Cj\[T^ can be optimized in N, 
so that the error term in the previous equation can be made exponentially small 
with respect to r. 

Such a result is very important and has many applications in the case where 
/ has some strong geometric properties, such as Hamiltonian or reversible struc- 
ture. In this situation, and under some compatibility conditions on the numerical 
method the modified vector field inherits the structure of /. For example 
if $r is symplectic and / Hamiltonian, then remains Hamiltonian. This has 
major consequences such as the preservation of a modified Hamiltonian over very 
long time (of order t~^) for the numerical solution, from which we can deduce 
long time stability results, existence of numerical invariant tori in the integrable 
case, etc... 

In the Hamiltonian case, this idea goes back to Moser [17], but was applied 
later to symplectic integrator by Benettin &; Giorgilli [3], Hairer &; Lubich [7] and 
Reich [18]. Such results now form the core of the modern geometric numerical 
integration theory for which we refer to the classical textbooks [8] and [14]. 
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More recently, these ideas have been extended in some situations to Hamilto- 
nian PDEs: First in the hnear case [4], and then in the semilinear case (nonhnear 
Schrodinger or wave equations), see [6, 5]. 

As far as stochastic differential equations (SDEs) are concerned, this approach 
has not been developed very much so far. Let us recall that given a SDE in M*^ 
of the form 

dX = f{X)dt + a{X)dW 

discretized by a numerical scheme - such as the Euler scheme for instance - with 
time step r providing a discrete sequence (Xp)pg then the error can be mea- 
sured in the strong or weak sense. Strong error means that Xp is a pathwise 
approximation of X^pr), and it is well known that the Euler scheme has strong 
order 1/2. Under standard assumptions on / and a, we have 

e( sup ||Xp -X(tp) in < ci (A;, T)r^/^ A; > 1, T > 0. 

\p=0,...,[T/r] J 

In this work, we consider another error which is often more important. We 
investigate the weak error which concerns the law of the solution. The Euler 
scheme has weak order 1. Under suitable smoothness assumptions on /, a and 
(/? : M'' ^ M (see for instance [10, 12, 20]): 

\E{^{Xp))-E{ip{X{tpm<C2{ip,T)T, p = 0,...,[r/r], T > 0. 

An attempt has been made by Shardlow [23] to extend the backward error analysis 
to this context. He has shown that the construction of a modified SDE associated 
with the Euler scheme can be performed, but only at the first step, ie for N = 2, 
and only for additive noise, ie when cr(X) does not depend on X. In this case, 
he is able to write down a modified SDE: 

dX = f{X)dt + a{X)dW 

such that 

E((^(Xp))-E((^(X(tp)))| <C3(^,r)r2, p = 0,...,[T/t], T > 0. 

He explains that for multiplicative noise or higher order, there are too many 
conditions to be satisfied by the coefficients of the modified equations. 

In this paper we take another approach, and build a modified equation not 
at the level of the SDE, but at the level of the generator associated with the 
process solution of the SDE. It is well known that given if : — )• R and denoting 
by X{t,x) the solution of the SDE satisfying X{0) = x, the function u{t,x) = 
E{cp{X{t,x))) satisfies the Kolmogorov equation 

dtu{t,x) = L{x,dx)u{t,x), 
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where L is the order 2 Kolmogorov operator associated with the SDE. 

In the case of the Euler method apphed to a SDE, we show that with the 
numerical solution, we can associate modified Kolmogorov operator of the form 

L{T,x,dx) = L{x,dx) +TLi{x,dx) + TL2{x,dx) H 

where L^, i > 1 are some modified operator of order 2i + 2. Again, the series 
does not converge but truncated series: 

,dx) = L{x,dx) + TLi{x,dx) H \-TLiy{x,dx) 

are considered. 

Note that in contrast with the classical case, we do not have a modified SDE 
and cannot straightforwardly define a solution to the modified equation 

dtv''{t,x) = L^''\T,x,d,)v''{t,x). 

However, in the case where the SDE is elliptic or hypoelliptic, we can build an 
approximated solution v^^^ such that 

E{ip{Xp))-v^^\pT,x)\<C4{ip,T,N)T^, p = 0,...,[T/t], T > 0. 

Furthermore, using the exponential convergence to equilibrium, we prove that 
in fact the constant C4 does not depend on T so that we have an approximation 
result valid on very long times. We also show that there exist a modified invariant 
measure for L(^)(r,j;,9^.). 

We can then use this weak backward error analysis to prove that the numer- 
ical solution Xp, p > 1 obtained by the Euler scheme is exponentially mixing up 
to some very small error, and for all times. This is typically a geometric numer- 
ical integration result in the sense that we prove the persistence of a qualitative 
property of the exact flow (exponential mixing) to the numerical approximation, 
over long times. 

Note that v^^^ is in fact constructed as a truncated series: v^^'' = X^^q '^^'^n 
and that vq = u is the solution of the Kolmogorov equation. Therefore, our 
result provides an expansion of the error as in [22] (see also [1, 2]). However, the 
expansion is different here. 

Error estimates on long times for elliptic and hypoelliptic SDEs have already 
been proved. In [15, 19, 20, 21], it is shown that for a sufficiently small time step 
the Euler scheme defines an ergodic process and that the invariant measure of 
the Euler scheme is close to the invariant measure of the SDE. In [22], the first 
term of an expansion of the invariant measure of the Euler scheme with respect 
to r is also given. In our work, we provide the expansion at any order. 

We emphasize that in our result, there is no particular smallness assumption 
on the stepsize r used to define the numerical solution. In particular, the discrete 
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process is not supposed to have a unique invariant measure, as in [15] or [19, 20, 
21]. 

This is also the case in the recent work [16]. There, it is shown that given an 
elhptic or hypoeUiptic SDE, the ergodic averages provided by the Euler scheme 
are asymptotically close to the average of the invariant measure of the SDE. 
Higher order schemes also considered. The main tool in [16] is the ellipticity or 
hypoellipticity of the Poisson equation, ie the equation L{x, dx)u = g. 

As in [16], we consider the case where the SDE is set on the torus T*^. This 
simplifies the presentation and the main ideas are not hidden by technical diffi- 
culties. In the same spirit, we only study the Euler scheme. In a forthcoming 
article, we will present more realistic applications of our method for SDEs set on 
R" with polynomial growth coefficients under suitable assumptions, and for more 
general schemes. As an example, we will treat the Langevin equation as in [21]. 



2 Preliminaries 

We consider the stochastic differential equation 

dX{t) = f{X{t))dt + a{X{t))dW, 

where the unknown X = (X*)j=i...^rf lives in the d- dimensional torus T*^. Also, 
/ = (/*)f=i ^iid a = {(T\{x))i=i^^^^^ii,i=i,...,m are smooth vector fields periodic in 
X G T'^'. The process {W^ {t) ^ . . . {t)) is a m-dimensional standard Wiener 
process over a probability space [Vi, J-", F) endowed with a ffitration {J-t)t>o- Using 
these notations, we rewrite the equation as: 

m. 

dx\t) = fix{t))dt + J2 <y}{x{t))dw'{t), i = i,...,d. 

In all the paper, smooth functions means C°° functions. Given a smooth 
function tp defined on T*^, we denote by HV'll^fc its norm in C*^(T'^,M). We also 
denote by ||'i/'||oo = sup^.g^pd For a multiindex k = {ki, . . . , kd) € N'^, we 

set \k\ = ki + ■ ■ ■ + kd and 

dxi^ ■ ■ ■ dx/ 

Therefore 

UWck ■= sup \dk'il^{x)\. 

3=Ui,---Jd) 
\j\<k 
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We also define the semi-norm 

IV'Icfc ■= sup |5fc^(x)|. 
j=(ii,---,id) 
i<|j|<fc 

In the fohowing, we assume that / and a are smooth and since we are work- 
ing with a stochastic difi^erential equation on the torus, standard theorems give 
existence and uniqueness of a solution for any initial data X(0) = x € T"^. We 
denote this solution by X{t, x), t >0. Also, since we chose to work on the torus, 
we do not have any problem of possible unbounded moments and this solution 
has clearly all moments finite. 

We denote by L{x,dx) the Kolmogorov generator associated with the stochas- 
tic equation: 

L{x,dx)v{x) = f{x)div{x) + a'^x)dijv{x), 
where we use the summation convention for repeated indices and di = dx^ and 

-. m 

«^'(^) := i^Y.^\{x)al{x). 

1=1 

It is well known that the Kolmogorov equation: 

^ = L{x, dx)u, X G T^, t > 0, u{0, x) = ip{x), x £ T'^, (2.1) 

with periodic boundary conditions has a unique solution for a smooth function ip 
and that for all x eJ'^: 

u{t,x) = E{ip{X{t,x))). 

Moreover, this solution is smooth. In the following, we write: u{t) = Pfcp so 
that {Pt)t>o is the transition semigroup associated with the Markov process 
(X(t,x))j>o xgT'*- Note that we use the standard identification u{t) = u(t, •). 

We wish to investigate the approximation properties of the Euler scheme for 
long times. We need assumptions on the long time behavior of the law of the 
solutions of (2.1), ie of the law of the Markov process. We assume the following 
mixing properties: 

[HI] There exists a C~(T'^,M) function p > such that 

L{x,dx)*p{x) = and / p{x)dx = 1. (2.2) 
In other words, the measure p{x)dx is invariant by X{t,x). 
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[H2] Let g € C°°(T"',R), and assume that J^.^ g{x)dx = 0. Then there exists a 
unique function //(x) G C°°(T'^,M) such that 

L{x,dxy fJ.{x) = g{x), and / fi{x)p{x)dx = 0. (2-3) 

[H3] Let u{t,x) be the solution of (2.1). Assume that fjj_ ip{x)p{x)dx = 0, then 
there exists a constant A and, for each /c G N, a polynomial Pk{t), such that 
if ip{x) € C°°(T'^,R) we have the estimates 

Vt>0, \\uit,-)\\ck<Pk{t)e-^'\Mck- (2.4) 

These hypothesis are usually satisfied under elliptic or hypoelliptic assump- 
tions on the operator L{x,dx)- The reader may find in [13] conditions to ensure 
[HI]. We also refer to [1, 16] for a general definition of hypoelliticity and applica- 
tions to numerical schemes. Combining kernel estimates for hypoelliptic diffusion 
([1, 11]) and exponential convergence to equilibrium, [H3] can be proved. Note 
that similar estimates are used in [15, 19, 20, 21], where specific examples are con- 
sidered. Finally, we mention that these hypothesis can be proved to be fulfilled 
using partial differential equations techniques (see [9]). 

For a smooth function ?/^, we set 

(ip) = / 'ip{x)p{x)dx. 

Note that by [H3], we have for any solution of (2.1) and A; € N 

Vt>0, \\u{t,-)-{ip)\\^,<p,{t)e~'^'y-{^)\\^,. (2.5) 

Now for a small time step r > and x S T"^, we consider the Euler method 
defined, for i = 1, . . . , by Xq = x and the formula 

= K + rfiXn) + al{X^){W'{{n + 1)t) - W^nr)), (2.6) 

for n > 0. Our main result can be stated as follows: 

Theorem 2.1 Let N and tq > be fixed. Then there exists a modified smooth 
density 

p^{x) = p{x) + Tpi{x) H h t'^ pn{x) 

such that Jijpd p^ {x)dx = 1, a constant Cn and a polynomial P/v(t) such that the 
following holds: For all smooth function function (p{x) on T'^, we have 

Vp G N, m^iXp) - j^^ (^d/i^ll^ < (P7v(ip)e-^*- + CjvT^) ||</^||c8iv+2 , (2.7) 

where for all p, tp = pr and dp^{x) = p^{x)dx. 
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This result can be viewed as a discrete version of (2.5). Note that it implies 
that all the invariant measure of the numerical process Xp are close to d//^ up 
to a very small error term C^t^ . 

Using this result, we can also recover the weak convergence result 

for some constant C{(p) depending on ip, and where we set dp{x) = p{x)dx. 
This can be compared with [15, 16, 19, 20]. As in [16], the only assumption 
made on r is that t < tq where tq is any fixed number. The influence of tq is 
only reflected in the constants in the right-hand side - we can for example take 
To = 1. In particular, we do not assume that Xp has a unique invariant measure 
- something that would be guaranteed only if r is small enough. We also recover 
an expansion of the invariant measure as in [22]. 

In the next sections, the constants appearing in the estimate depend in general 
on bounds on derivatives of / and g defining the SDE. They will also depend in 
general on tq and N, but not on ip. 

3 Asymptotic expansion of the weak error 

We have the formal expansion for small t: 

M(t, x) = ip{x) + tL(x, dx)^{x) + —L{x, dxf(pix) H h —:L{x, dxTip{x) ^ 

2 n\ 

This is just obtained by Taylor expansion in time since by (2.1) ■^u{t,x) = 
L{x,dx)^u{t,x) and in particular ■^u{0,x) = L{x,dx)"'(p{x). 

Since the solution u of the Kolmogorov equation is smooth and has its deriva- 
tives bounded in terms of the initial data ip, the above formal expansion can be 
justified and we have the following proposition whose proof is easy and left to the 
reader: 

Proposition 3.1 Assume that (p G C°°(T'^,M), and let tq > 0. Then for all N, 
there exist a constant Cn such that for all r < tq, 

^ n 

- E ^^(^,5x)Xx)||^ < C;vr^+Vlc2iv+2 • (3.1) 

n=0 

With the Euler scheme defined in (2.6), we associate the continuous process 

Xi{t) = Xl + {t- nT)f\Xn) + a\{Xn){W\t) - W^nr)), t € [nr, {n + 1)t], 

(3.2) 
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and X(nr) = X„. We thus have = Xx{tn+i)- The process (3.2) satisfies 

the equation 

dXiit) = fiXn)dt + a}{Xn)dW\t), t £ [nr, {n + l)r]. (3.3) 

Clearly, (X„) defines a discrete in time homogeneous Markov process but X^ is 
not Markov. 

In this work, we are only interested in the distributions of the solutions and 
of their approximation. We now examine in detail the first time step and its 
approximation properties in terms of the law. By Markov property, it is sufficient 
to then obtain information at all steps. Next result gives an expansion similar to 
Proposition 3.1 for the Euler process. 

Theorem 3.2 Then for all n > 1, there exist operators An{x,dx) of order 2n, 
such that for all N >1, there exist a constant Cm satisfying 

N 

||Ev9(X,(t)) - J^rM„(x,a,)(^(x)||^ < C;vr^+Vlc2iv+2 , (3-4) 

n=0 

for all ip e C°°(T'=', M) and x S T'^. 

Proof. Using (3.3) and the Ito formula, we get for t < t, 

dip{Xx{T)) = L{x,dMXx{t)) + al{x)^^ip{Xx{t))dW\x), 
or equivalently, 

ifiXxit)) = ifix) + f L{x,dMXx{s))ds+ [' ai{x)dMXAs))dW\s). (3.5) 
Jo Jo 

Note that the last term is a martingale. We define the operator 

Ro,i{x,d.x) = (T\{x)di. 

We have for all s and all x, 

L{x,dMXx{s)) = {f\x)diip){Xx{s)) + {a'^{x)dijip){Xx{s)). 
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Hence applying (3.5) to Sjy? and dij<p, we obtain 
L{x,dx)ip{Xx{s)) = L{x,dx)'f{x) 



+ / r(x)/^'(x)a,,(/.(X,.(a)) + r(x)a"^(x)a„,,^(X,.(a))d(T 
Jo 







+ / f\x)al{x)d,j^{X,{a))dW'{a) 



r a'^ ix)ai ix)dinM^^ {a))dW\a) . 
Jo 
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We set Aq = I, Ai = L and plugg this in (3.5) to obtain 

ip{X:^{t)) = ip{x) +tAi{x,dx)^p{x) + / A2{x,dx)ip{Xx{cr))dads 

Jo Jo 

+ [ [ Ri,e{x,dMXAc^))<iW\a)ds + [ Ro,e{x,dMXA<^))dW'{a). 
Jo Jo Jo 

where 

A2ix, d,) = f{x)P{x)d^j + r{x)a^\x)dnu + a'\x)r{x)d^,n + a^^ (x)a"^(x)%„, 

is an operator of order 4, and 

Ri,i{x,d^) = f{x)al{x)dij + a''^{x)al{x)dinj 

are operators of order 3. Taking the expectation so that the last two term disap- 
pear, we easily deduce the result for = 1. 

Let us now prove recursively that there exist operators An{x, dx) of order 2n and 
Rn/{x, dx) of order 2n + 1, such that 

AT 

ifiXxit)) = ip{x) +tL{x,dxMx) + J2t''Mx^9Mx) 

n=2 

+ •••/ AN+iix,dx)ipiXxisN+i))dsi- ■ -dsN+i 
Jo Jo 

+ ■■•/ Rn,eix,dx)ip{Xx{sn+i))dsi- ■ -dSndW^Sn+l). 

n=l-^0 Jo 

(3.6) 

Note the expectation of the last term vanishes so that (3.6) easily implies (3.4). 
To prove (3.6), assume that ^at+i and Rj\f/ are known, and let us decompose 
An+1 as AN+iix,dx) = Al^_^_^{x)dj, where j = (ji,...,jm) are multiindices 
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(with the summation convention) and smooth functions of x. Such decom- 

pistion is easy to write for = 1 or 2. We apply (3.5) to djip{Xx{sN^i)) for a 
given muti-index j, and obtain 

f-SN+l 

+ / A'j^,^{x)P{x)dndjip{X^isN+2))dSN+2 

Jo 



f-SN+l 

+ / A]^,^{x)a"'\x)dn£djip{X^{sN+2))dSN+2 

Jo 

+ / A'j^^^{x)anx)dndjip{X^{sN+2))dW\sN+2). 
Jo 



We thus choose 

An+2{x) = A%^^{x)r{x)dndj + A%^^{x)a'^{x)d,jdj. 

and 

Rn+i,i{x) = A^j^^^{x)a1{x)dndj, 
and obtain (3.6) with + 1 replaced by + 2. ■ 

4 Modified generator 
4.1 Formal series analysis 

Let us now consider r as fixed. We want to construct a formal series 

L(t; X, d^) = L{x, d^) + tLi{x, d^) + ■ ■ ■ r"L„(x, d^) + ■ ■ ■ (4.1) 

with operator coefficients Ln{x,dx) smooth on T'^, and such that formally the 
solution v{t, x) at time t = r of the equation 

dtv{t, x) = L{t; X, dx)v{t, x), u(0, x) = ip{x) 

coincides in the sense of asymptotic expansion with the approximation of the 
transition semigroup E,{lp{Xx{t))) studied in the previous section. In other words, 
we want to have the equality in the sense of asymptotic expansion in powers of r 



exp(rL(r; x, dx))ip{x) = ip{x) + ^ r'^A„(x, dx)ip{ 
where the operators An{x,dx) are defined in Theorem 3.2. 
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Formally, this equation can be written 

exp(rL(r; x, d^)) - Id = rl(r) (4.2) 

where I(r) = Y.n>i ^""'^n- 
We have 

exp(TL(r; x, S^^)) - Id = tL{t; x, dx)[ - — — — L(r; x, 5^,; 

n>0 ^ ' 

Note that the (formal) inverse of the series is given by 

n>0 ^ n>0 

where the are the Bernoulli numbers: see for instance [8, 5] and [4] for a 
similar analysis involving operators. Hence equations (4.1), (4.2) are equivalent 
in the sense of formal series to 



l>0 

n 



n>0 y £=1 niH \-ni+nij^i=n—t I 

(4.3) 

Identifying the right hand sides of (4.1) and (4.3), we get the following recursion 
formula 

Ln = An+l + Y,^ Yl Ln,---Ln,An,^,+l. (4.4) 

£=1 niH \-ni+nm_i=n—£ 

Each of the terms of the above sum is an operator of order 2n + 2 with smooth 
coefficients and therefore L„ is also an operator of order 2n + 2 with smooth 
coefficients. 

Note that (4.2) gives immediately the inverse relation of this formal series 
equation: 

n ^ 

^- = ^£1 Yl Ln,---Ln,. (4.5) 

(.=1 ni+-+ni=n—i 

Moreover, we have clearly 

Ln{x,d^)l = 0, 

where 1 denote the constant function equal to 1. 
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4.2 Approximate solution of the modified flow 

For a given A^, we have constructed in the previous section an operator 

N 



L^^\T-x,d^) = L{x,d.,) + Y,T^Ln{x,d.,). (4.6) 

n.= l 

In order to perform weak backward error analysis and estimate recursively the 
modified invariant law of the numerical process, we should be able to define a 
solution v^{t,x) of the modified flow 

dtv''{t,x) = L^''\T;x,d,)v''{t,x), z;^(0,x) =v?(x). (4.7) 

However, in our situation we do not know whether this equation has a solution. 
This is in contrast with standard backward error analysis where the modified flow 
can always be defined. 

The goal of the next proposition is to give a proper definition of the modified 
flow (4.7). 

Theorem 4.1 Let if be a smooth functions on T^. For all n € N, there exist 
smooth functions V£{t,x), defined for all times t >0, and such that for all t > 
and n € N, 

n 

dtVnit, x) - LVn{t, x) = ^ L^f^.f (t, x), (4.8) 

e=i 

with initial conditions vq{0, x) = (p{x) and u„(0, x) = for n > 1. For all N > 0, 
setting 

N 

v^^Ht,x)=J2r''vn{t,x), (4.9) 

k=0 

then the following holds: 

(i) There exists a constant Cm such that for all time t > 0, and all r > 0, 

||Ez;(^)(t,l,(r))-t;(^)(t + r,x)||^ <CjvT^+i sup \vn{t + s)\ . (4.10) 

sG{0,t) 
n=0,...,N 

(ii) For fixed tq > 0, there exists a constant Cn such that for all r < tq, 

||E<^(X.(r)) -t;(^)(r,x)||^ < Cmt^'^'M^,,^, . 

Proof. For n = 0, the equation (4.8) implies VQ{t,x) = u{t,x), the solution of 
(2.1). Let n € N and assume that Vj{t,x) are constructed for j = 1 . . . ,n — 1. 
Let 



Fn{t,x) ■.= Y,LiVn-i{t,x) (4.11) 
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be the right-hand side in (4.8). Then f„ is uniquely defined and given by the 
formula ^ 

Vnit)= [ Pt^sFnis)ds, t > 0. (4.12) 

Jo 

By [H3], it is not difficult to check that v^, n € N are smooth and that for all 
t > 0, 

\\vnmc,<c{t)M^,^,^, ken, (4.13) 

where the constant C{t) depends on t, k, n and on the coefficient of the equation. 
Clearly, this type of estimate can be improved in the elliptic case. This proves 
the first part of the Theorem. 

To prove (i), we consider a fixed time t, and define the functions Wn{x,s) := 
Vnit + s) for s > and n € N. By definition, these functions satisfy the relation 

n 

dsWn{s,x) - LWn{s,x) = LgWn-iis , x) , Wn(0,x) = Vn{t,x). 

e=i 

Let us consider the successive time derivatives of the functions Wn{s, x). We have 
using the definition of Wn, for all s > 0, 

n n 

d^Wn{s, x) = ^ LidsWn-e{s, a;) = ^ ^ Li^L^^Wn-kis, x), 

and we see by induction that for all m > 1 and s > 

d'^Wn{s,x)= ^ Le^---Le^we^^,{s,x). (4.14) 

\-im+i=n 

Using the fact that the operators Le are of order i + 2 with no terms of order 
zero, we see that there exists a constant C depending on n and m, such that 

\\d^Wnis,x)\\^ < C ^sup kfc(s)|^2,+2„ . 
Now let us consider the Taylor expansion of Wn{T), for t < tq and n = 0, . . . , N , 

wnir,x) = J^L^ar(0,x)+ / -af— +iu;„(.,x)ds 

^^ml Jo {N-m)l 

N-n 

= y] — T Y] Li^---Li^Wi {0,x)+RN,niT,x). 

^ — ' ml ^ — ' 

m=0 ^iH \-e,n+i=n 

Using the bounds on the time derivatives of Wn{s, x), we obtain that for all r > 
and all n = 0, . . . , A^, 

||^7V,n(-r,x)||^ < Cr^-™+^ sup k„(s, x)!^^^^^ , 

sG(0,r) 
n=0,...,N 
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for some constant depending on N, m. After summation in n, and using the 
expression (4.5) of the operators An and the definition of Wn, we get 

N n 

V^^^t + T,x) = X]^" AniVn-m{t,x) + RN{t,T,x) 

where 

\\RN{t,T,x)\\^ < CNr^~^^ sup |z;n(t + S,x)|^4jv+2 ■ 

se(0,r) 
n=0,...,N 

To conclude, we use (3.4) applied to (f = v^^\t,x), and we easily verify that 
E?;W(t,X^(t)) satisfy the same asymptotic expansion. 

The second estimate (ii) is then a consequence of (i) with t = and (4.13). ■ 

Note that in the previous theorem, we have constructed a function v^^\t,x) 
which is an approximate solution of (4.7). More precisely, we can easily show 
that we have for all time t >0, 

dtv^^\t,x) = L^^\T;x,d^)v^^^{t,x) + R^^\t,x), v^^\o,x) = ^{x), 

where 

R^''\t,x) = - Yl r'^^''L,,v,,{t,x) 
ei,e2=o,...,N 

ll+t2>N 

is of order 0{t^+^). 

5 Asymptotic expansion of the invariant mea- 
sure and long time behavior 

We now analyze the long time behavior of the solution of the modified equation 
(4.7). In the following, for a given operator B{x,dx), we denote by B{x,dx)* its 
adjoint with respect to the product. We start by an asymptotic expansion of 
a formal invariant measure for the numerical scheme. 

Proposition 5.1 Let (Ln)n>o be the collection of operators defined recursively by 
(4.4). There exists a collection of functions (/x„(x))„>o such that Ho{x) = p{x), 
Jjd Mr!,(3;)dx = for n > 1, and for all n > 0, 

n 

Lol^n = -'^iLi)* fin-e- (5.1) 
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Let N > be fixed and LW(r;x, dx) the operator defined by (4.6). Then the 
function 



N 



(r;x) = + ^r"/i„(x) G C°°(T^] 



n=l 



satisfies 

[ /i(^)(r;a;)dx = 1, 

and 

LW(r;x,a.)V^)(r;x) = G(^)(r;x), 
with for all k and all r € [0, tq], 

||GW(r;x)|L, < C^fcT^+i and [ G^^) (r; x)dx = 0, 

where Cjsi,k depends on N and k. 

Proof. Assume that fiQ = p and fij{x) are known, for j = 0, ... ,n — 1 with 
n > 1. Consider the equation (5.1) given by 

n 

-^oMn = ~ ''^{Li)* Hn-i = Gn- 
1=1 

Note that the right-hand side Gn{x) is a smooth function satisfying 

/ G„(3;)dx = - / (L£)*/i„_£dx = - / /z^.^L^ldx = 0, 

where 1 denotes the constant function equal to 1, which as already seen is in the 
Kernel of all the L^. 

Using Hypothesis [H2], we easily obtain the existence of a C°° function /u„ satis- 
fying (5.1) and Jj^ lJ^n{x)dx = 0. This shows the first part of the Proposition. 
We then write 

27V 

LW(T;x,a,)V(^)(T;x) = J^r" ^ (L^JV^, 

n=0 ^i+£2=n 
2N 

n=N+l ei+£2=n 

e,<N 

=■■ G(^)(r;x), 

and we easily verify that G^^^ satisfies the hypothesis of the Proposition, owing 
to the fact that L/1 = for aU i > 0. ■ 
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Proposition 5.2 For all n and k there exists a polynomial Pk,n{t) such that for 
all t > 0, 

\\Vn{t,x) - / ifdHnhk < ^'fc,„(t)e"^*||v3 - {^)hk+4„ ■ (5-2) 

Proof. Using the fact that //q = p and vq = u, wee see that estimate (5.2) 
is satisfied for n = (see Equation (2.5)). Let n > 1 and assume that vj, j = 
0, . . . , n - 1 satisfy for A; € N, i > 0: 

\\vj{t,x) - / ip{x)nj{x)dx\\ < Pkj{t)e"^^\\ip{x) - , 

for some polynomial Pkj- 
Let us set 

n ,. 

Cn(t) = / Vn-m{t,x)lJim{x)dx. 

Td 

m=0 ^ 

We claim that Cn{t) does not depend on time. Indeed: 

n „ n „ 

y^dt I Vn-m{t,x)^lm{x)dx = ^^Ot Urn{t, x)fln-mix)dx 

Td Td 

m=0 ^ m=0 '' ^ 

n m „ 

= / Lm-eve{t,x)fin-m{x)dx 

m=Oi=0 J'^' 
n n „ 

= VV / Vi{t,x)L*^_f^lin-m{x)dx 

n „ n—l 
= / Vl{t,x)y^L*^l_Ln-ll-m{x)dx = {), 

Td 

1=0 •' ^ m=0 

by definition of the coefficients see (5.1). Note that, thanks to the smoothness 
properties of all the functions, the computation above is easily justified. 
We deduce: 

/ dtVn{t,x)p{x)dx = -S^ \ dtVn-~m{t,x)fIjnix)dx. (5.3) 
Jjd Jjd 

Next, we compute the average of F„. By (4.8), (4.11) and (5.3), we have 



{Fn{t)) = I Fn{t,x)p{x)dx = / dtVnit,x)p{x)dx - / Lvnit,x)p{x)dx 

Ifd Jjd 

dtVn{t,x)p{x)dx 

Jd 



n „ 

> , / dtv 

n—m 

{t,x)fim{x)dx. 

Td 



n 



We rewrite (4.12) as follows 

Vn{t,x)= [ {Fn{s))ds+ [ (F„ (s, x) - (F„ (s)) )ds. 

Jo Jo 

Using the previous expression obtained for (F„,(s)) and recalling the initial data 
for Vn, we deduce that 



Vn{t,x) = - S2 / Vn-m{t,x)Hrn{x)dx + / 
m=l ■^T'* 



t 



ip{x)fj,n{x)dx 
+ / Pt-s{Fn{s,x)-{Fn{s)))ds. 

Jo 

Then, using Hm{x)dx = 0, m G N, we get 

Vn{t,x)- ip{x) Hn{x)dx = y2 / {Vn-m{t,x)- (f{x) fln~m{x)dx) Hrnix)dx 

+ / Pt-s{Fn{s,x)-{Fn{s)))ds. 



Jo 

Note that, since L^, £ G N is a differential operator of order 2i + 2 with smooth 
coefficients and containing no zero order terms, we have 

n-l 
£=0 

Then by [H3] 



n— 1 „ 

||F„(s) - (i^„(s))||^fc < Cfc,^^ - / Vi{s,x)fj,i{x)da 

o—c\ J"^"^ 



\ Qk + 2(n-e)+2 ■ 



\Vn{t,x)- / Lp{x)lln{x)dx\\^< 

Jjd 



\flni{x)\dx 

jd 



\\Vn-m{t,x) - / ip{x)lJ,n-m{x)dx\\ / 

+ r - s)e-^(*-^) ||F„(S, x) - 11^, ds, 



and using the recursion assumption 

„ n 
\\Vnit,x) - / ip{x)fJ,nix)dx\\^^ < ^ (t)e"*^||9.(x)-(^)|L 

m=l 

n-l ..t 

+ 5Z / P'^^* ~ 5)^'fc+2(n-£)+2/(s)e"^*ds||(/j(x) - (99)11^^ 
l—n Jo 



e=o 

The conclusion follows 
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We give now our main result concerning the long time behavior of the numer- 
ical solution: 



Theorem 5.3 Let tq and N be fixed. Then there exists o,nd a polynomial 
PN{t) such that the following holds: Let Xp be the discrete process defined by 
(2.6); then we have for p >0, r > and smooth function (p{x) 

VpGN, \\E^{Xp)-v^''\tp,x)\\^<CNT^M\csN+2, (5.4) 
where for all p, tp = pr. Moreover, we have 

Vp G N, \\E^{Xp) - I '^d^^'^^IL < (PN{tp)e-^'^ + Cnt 
where dfi^^\x) = fi^^\x)dx. 



'P\\(2SN+2 5 

(5.5) 



Proof. For all p, with tj = jr, we have 

EyD{Xp)-v^^^{tp,x) = Ev^^\0,Xp) -v^^\tp,x) 

p-i 

= E^E^^--(z;(^)(t,-,Xp_,-) -^(^)(t,+i,Xp_,-„i)^ 

Here we have used the notation E^p-^-'^ for the conditional expectation with 
respect to the filtration generated by By the Markov property of the 

Euler process at times tj: 

^yiN)(^t^^Xp_j) - t;(^)(t,+i,Xp_,-_i)) 

= E^— ^ (^(^) (t, , (r)) - ^;(^) X,_,_i)) . 

Using (4.10) with t = tj, and Proposition 5.2, we deduce that 

p-i 

\\E^{Xp)-V^''\tp,x)\\^ < C^T^+l^ sup \Vn{tj + S,x)\^^^^ 

j^Q s£(0,t) 



n=0,...,N 
p-1 



j=0 

for some constant Cn and polynomial (5iv(i)- We have used: \vn{tj + s, x)\^j^_^_2 ~ 
{''^nitj + s,x) — jjd '^dfin\^j^_^2' Conclude by using the fact that for a fixed 
constant 7 > 0, we have 

ye-7.-< 1 <^ 

- 1 _ - r 

j=0 
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where the constant C depends on 7 and tq. This shows (5.4). The second estimate 
is a consequence of Proposition 5.2. ■ 
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